********************************************************************************
*   Title: draft_kinky.do
* Purpose: Produces final graphs and tables for main draft, excluding appendix.
*          See the main() program at the bottom for organization.
* Authors: Qiping Xu, Eric Zwick
*     RAs: Tom (Tianfang) Cui, Laurence O'Brien, Harleen Kaur,
* 		   Francesco Ruggieri, Caleb Wroblewski, Ramiro Rossi
*    Date: 2024-05-02
********************************************************************************

******************************************************************
* Figure 1: Time Series of Fiscal Q4 Investment Spikes (1984-2016)

capture program drop fig_sawtooth /* %< */
program define fig_sawtooth
	
	*plot a quarterly CAPX normalized by each firm, and then average across firms
    *to generate a time-series sawtooth pattern 
	
	* Prep data (see programs_kinky.do)
	
	qvar_prep capx, auxvars(at cal_month sale4_3 oibdpy3 oibdpq4 cfy3 cfq4 ///
        capxq4 capxqf11)

	* Figure 1a	
	
	preserve
		qvar_plot, varlab(CAPEX) output(sawtoothpattern_p50.pdf) ///
		ysc(75 125) xsc(96 219) xlab(99 (40) 219) xsize(7)
	restore	
	
	* Figure 1c
	
	preserve
		drop if cal_month==12
		qvar_plot, varlab(CAPEX) output(stpnondec.pdf) ///
		ysc(70 125) xsc(96 219) xlabel(99 (40) 219)  
	restore
	
	* Figure 1d
	
	preserve
		gen avecf = cfy3 / 3
		drop if cfq4 > avecf
		qvar_plot, varlab(CAPEX) output(stpcf.pdf) ///
		ysc(70 125) xsc(96 219) xlabel(99 (40) 219)  
	restore	
	
end /* %> */

capture program drop fig_sawtooth_eps /* %< */
program define fig_sawtooth_eps
	
	*plot a quarterly CAPX normalized by each firm, and then average across firms
    *to generate a time-series sawtooth pattern 
	
	* Prep data (see programs_kinky.do)
	
	qvar_prep capx, auxvars(at cal_month sale4_3 oibdpy3 oibdpq4 cfy3 cfq4 ///
        capxq4 capxqf11)

	* Figure 1a	
	
	preserve
		qvar_plot, varlab(CAPEX) output(sawtoothpattern_p50.pdf) ///
		ysc(75 125) xsc(96 219) xlab(99 (40) 219) xsize(7)
	restore	
	
	* Figure 1c
	
	preserve
		drop if cal_month==12
		qvar_plot, varlab(CAPEX) output(stpnondec.pdf) ///
		ysc(70 125) xsc(96 219) xlabel(99 (40) 219)  
	restore
	
	* Figure 1d
	
	preserve
		gen avecf = cfy3 / 3
		drop if cfq4 > avecf
		qvar_plot, varlab(CAPEX) output(stpcf.pdf) ///
		ysc(70 125) xsc(96 219) xlabel(99 (40) 219)  
	restore	
	
end /* %> */

capture program drop fig_quarterly_capx /* %< */
program define fig_quarterly_capx

	* Figure 1b
	
	drop if capxq==.|capxq<0
	drop if fyear<1984|fyear>2016
	winsor2 capxq, cut(1 99) replace 
	keep gvkey fyear fyearquarter fquarter capxq  
	tsset gvkey fyearquarter
	rename fyearquarter fqtr
	format fqtr %tq
	preserve
	collapse (median) capxq, by(fqtr)
	forv j=99(4)239 {
		local spikes `spikes', `j'
	}
	
	twoway (line capxq fqtr) (scatter capxq fqtr if inlist(fqtr `spikes'), mcolor(red)), ///
		ytitle("Median Quarterly CAPEX (\$ Mils)") xtitle("Fiscal Quarter") ///
		xsc(range(96 219)) xlabel(99 (40) 219) ysc(range(1 6)) ylabel(1 (1) 6) ///
		graphregion(color(white) lwidth(large)) bgcolor(white) legend(off) 
	graph export "quarterlycapxsample.pdf", as(pdf)  replace

end /* %> */

capture program drop fig_fyrchange /* %< */
program define fig_fyrchange
	
	* Figure 1e
		
	bysort gvkey: egen maxdif=max(dif)
	keep if maxdif==18

	drop if capxq<0
	gen capxq2ave=capxq/(capx/4)
	bysort gvkey fyear: egen sumratio=sum(capxq2ave)
	keep if sumratio==4
	drop sumratio
	
	gen start=yearquarter if change_ind==1&fquarter==4
	bysort gvkey: egen maxstart=max(start)
	drop quarter 
	gen quarter=yearquarter-maxstart
	keep gvkey fyear fyearquarter yearquarter fyrend_month quarter change_ind capxq2ave
	drop if quarter<-13 |quarter>4
	bysort gvkey: gen count=_N
	drop if count<12
	drop count

	collapse (mean) capxq2ave, by(quarter)
	sort quarter
	replace capxq2ave=capxq2ave*100

	twoway (bar capxq2ave quarter, color(blue)) ///
		(bar capxq2ave quarter if quarter==-10|quarter==-6|quarter==-2|quarter==2|quarter==6, ///
		bstyle(legend)) (bar capxq2ave quarter if quarter==0|quarter==4, ///
		bstyle(sunflowerdb)), xlabel(-12 (4) 4) ytitle(" CAPEXQ/Ave(CAPEXQ) %") ///
		xtitle("Quarter") yscale(range (70  130)) ylabel(70 100 130) ///
		graphregion(color(white) lwidth(large)) bgcolor(white) legend(off)
	graph export "fyrchange.pdf", as(pdf) replace
	
end /* %> */

******************************************************************
* Figure 2: International Evidence of Fiscal Q4 Spikes (2004-2016)

capture program drop fig_sawtooth_intl /* %< */
program define fig_sawtooth_intl

    rename capxy4 capx
    qvar_prep capx, auxvars(fic)
    bys fic: qvar_plot, varlab(CAPEX) output(sawtooth_intl) ///
		ysc(60 140) xsc(176 219) xlab(#5) ylab(80(20)140)

end /* %> */

***********************************************
* Figure 3: Fiscal Q4 Spikes and Tax Incentives
* Note: These figures utilize internal data
	
capture program drop fig_spike_bandwidth /* %< */
program define fig_spike_bandwidth

    load_analysis_data taxable_rd
	
	* Figure 3a
    #delimit ;
    twoway (scatter capx4_3 ni_bins, symbol(O)),
		xline(0, lcolor(gs9)) $graphconfig
        yline(100, lcolor(red))
        xtitle("Taxable Income Before Depreciation Bin ($000s)") 
        ytitle("Median CAPX 4/3");
    #delimit cr
    graph export "taxable_rd.pdf", replace

    load_analysis_data taxable_nolstock

	* Figure 3b
    #delimit ;
    twoway (scatter capx4_3 nol_group, symbol(O)) (lfit capx4_3 nol_group),
		$graphconfig legend(off)
        xtitle("Bin of Loss Carryforward to Net Income") 
        ytitle("Median CAPX 4/3");
    #delimit cr
    graph export "taxable_rd_nolstock.pdf", replace

end /* %> */

*****************************************************************
* Figure 4: Cross-Sectional and Dynamic Determinants of Q4 Spikes

capture program drop fig_uncertainty /* %< */
program define fig_uncertainty 

	* Figures 4c and 4d 

	load_analysis_data capx_trim5_2020

    local winspct 0.01
    local salesfilter 300

    uncertainty_prep, winspct(`winspct')
    replace capx4_3 = capx4_3/100

    foreach v of varlist ebitda2atsd ebitda2atmean {
        xtile `v'_bin=`v', nq(100)
        xtile `v'_bin30 =`v', nq(30)
    }

    * Mean
    gen loss = ebitda2at <= 0
    by gvkey: egen loss_share = mean(loss)
    replace loss_share = 100 * loss_share

    * Variance
    gen ebitda2at_pos = ebitda2at * (ebitda2at > 0)
    gen ebitda2at_neg = ebitda2at * (ebitda2at <= 0)
    replace ebitda2at_pos = . if ebitda2at == .
    replace ebitda2at_neg = . if ebitda2at == .
    by gvkey: egen ebitda2atsd_all = sd(ebitda2at)
    by gvkey: egen ebitda2atsd_pos = sd(ebitda2at_pos)
    by gvkey: egen ebitda2atsd_neg = sd(ebitda2at_neg)
	
    gen test = ebitda2atsd_all^2 - ebitda2atsd_pos^2 - ebitda2atsd_neg^2
    sum test, d
    gen variance_negshare = ebitda2atsd_neg / ebitda2atsd_all
    replace variance_negshare = 100 * variance_negshare
    tempfile stats
    save `stats'

    use `stats', clear
    collapse (p50) capx4_3 (mean) loss_share (mean) ebitda2atmean_bin, ///
		by(ebitda2atmean_bin30)

    #delimit;
    twoway 
        (scatter capx4_3 ebitda2atmean_bin, mcolor(navy)) 
        (lfit capx4_3 ebitda2atmean_bin, lc(maroon)),
        ytitle("Median CAPEX 4/3")
        yscale(range(100 140)) ylabel(100 (10) 140) 
        xtitle("Bin of EBITDA Mean") $gpr legend(off);
    graph export "ebitda2atmean.pdf", replace;

    twoway
        (bar loss_share ebitda2atmean_bin, color(navy*0.6) barwidth(2)),
        ytitle("Loss Share (%)")
        yscale(range(0 100)) ylabel(0(50)100)
        xtitle("Bin of EBITDA Mean") $gpr legend(off)
        ysize(2) scale(2);
    graph export "ebitda2atmean_loss.pdf", replace;
    #delimit cr

    use `stats', clear
    collapse (p50) capx4_3 (mean) variance_negshare (mean) ebitda2atsd_bin, by(ebitda2atsd_bin30)

    #delimit;
    twoway 
        (scatter capx4_3 ebitda2atsd_bin, mcolor(navy)) 
        (lfit capx4_3 ebitda2atsd_bin, lc(maroon)),
        ytitle("Median CAPEX 4/3")
        yscale(range(100 140)) ylabel(100 (10) 140) 
        xtitle("Bin of EBITDA Volatility") $gpr legend(off);
    graph export "ebitda2atsd.pdf", replace;

    twoway
        (bar variance_negshare ebitda2atsd_bin, color(navy*0.6) barwidth(2)),
        ytitle("Variance from Loss (%)")
        yscale(range(0 100)) ylabel(0(50)100)
        xtitle("Bin of EBITDA Volatility") $gpr legend(off)
        ysize(2) scale(2);
    graph export "ebitda2atsd_loss.pdf", replace;
    #delimit cr

end/*%>*/

	* Figure 4c
capture program drop fig_autocorrelation /* %< */
program define fig_autocorrelation
	
	replace capx4_3 = capx4_3 / 100

	xtset gvkey fyear
    sort gvkey fyear

    gen spike = capx4_3 > 134

    forv i = 1/10 {
        xi: areg F`i'.spike spike i.fyear, absorb(gvkey) cluster(gvkey)
        gen beta`i' = _b[spike]
        gen se`i' = _se[spike]
    }
    keep if _n == 1
    reshape long beta se, i(gvkey) j(lag)

    set obs `=_N + 1'
    replace lag = 0 if lag == .
    replace beta = 0 if beta == .
    replace se = 0 if se == .
    sort lag

    gen beta_U = beta + 1.96*se
    gen beta_L = beta - 1.96*se
    
    #delimit
    twoway 
        (rcap beta_U beta_L lag)
        (connect beta lag, msym(o) mc(navy) lc(navy) lp(dash)), 
        yline(0, lc(red)) 
        xsca(range(-.5 10.5)) xlab(0(1)10)
        ysca(range(-.09 0.01)) ylab(-.08(.02)0)
        ytitle("P(Future Spike | Spike in Year 0)") 
		ylabel(, format(%04.2f))
        xtitle("Years from Spike") $gpr $graphconfig legend(off);
    #delimit cr
    graph export "prob_future_spikes.pdf", replace

end /* %> */


* Figure 4d: Cumulative Investment after Q4 Spikes.

capture program drop fig_cumlevel_qu /* %< */
program define fig_cumlevel_qu

	su eventtime
	drop if eventtime > 12 | eventtime < 1
    gen avepre1ylow = avepre1y - 1.96 * sdavepre1y 
    gen avepre1yup = avepre1y + 1.96 * sdavepre1y 
    su avepre1yup
	local ymax = r(max)
    su avepre1ylow
	local ymin = r(min)
	
	#delimit ;
	graph twoway 
		(connect avepre1y eventtime if eventtime > 0,
			msize(medium) mcolor(navy) xtitle("Quarter") yline(100, lc(gray))
			xline(4, lc(gray) lp("-")) ysize(3.5)
			ytitle("Cumulative Average CAPEX Relative to Baseline (%)") 
			legend(off) xlabel(1 2 3 4 5 9) ysc(range (`ymin' `ymax')) 
			ylabel(50(50)250) graphregion(color(white) lwidth(large)) bgcolor(white))    
		(connect avepre1y eventtime if eventtime < 2, msize(small) lp(dot) mcolor(navy)) 
		(line avepre1ylow eventtime, lpattern(dot) sort)
		(line avepre1yup eventtime, lpattern(dot) sort),
		name(ave4, replace);
	graph export "aveq4spikerbasepre1y.pdf", as(pdf) replace;
	#delimit cr

end /* %> */

**********************************************
* Figure 5: Bunching
* Note: These figures use internal data 

capture program drop fig_bunching /* %< */
program define fig_bunching

    local graphopts graphregion(color(white) lwidth(large)) bgcolor(white) 
    
    ****************************************************************************
    * 1. Volatility graphs
    ****************************************************************************
    load_analysis_data bunch_mean
    replace buncher = 100 * buncher
    sort early ebitda2atmean_bin
    by early: gen bunch_rel = 100 * buncher / buncher[1]

    #delimit ;
    twoway (scatter bunch_rel ebitda2atmean_bin if early == 1) 
        (lfit  bunch_rel ebitda2atmean_bin if early == 1, lc(red) lp(dash))
        (scatter  bunch_rel ebitda2atmean_bin if early == 0, ms(S) mc(maroon)) 
        (lfit bunch_rel ebitda2atmean_bin if early == 0, lc(red) lp(dash)),
        ytitle("P(Bunch at Threshold, Indexed in % to Bin 1)")
        /*yscale(range(100 140)) ylabel(100 (10) 140)*/
        xtitle("Bin of EBITDA Mean") `graphopts'
        legend(order(1 "Early" 3 "Late"));
    graph export "ebitda2at_mean_bunch.pdf", as(pdf)  replace;

    gen ebitda2atmean_bin2 = ebitda2atmean_bin - 0.5 if early == 1;
    replace ebitda2atmean_bin2 = ebitda2atmean_bin + 0.5 if early == 0;

    twoway
        (bar loss_share ebitda2atmean_bin2 if early == 1, color(navy) barwidth(1))
        (bar loss_share ebitda2atmean_bin2 if early == 0, color(maroon) barwidth(1)),
        ytitle("Loss Share (%)")
        yscale(range(0 100)) ylabel(0(50)100)
        xtitle("Bin of EBITDA Mean") $gpr legend(off)
        ysize(2) scale(2);
    graph export "ebitda2atmean_loss_bunch.pdf", replace;
    #delimit cr

    load_analysis_data bunch_sd
    replace buncher = 100 * buncher

    sort early ebitda2atsd_bin
    by early: gen bunch_rel = 100 * buncher / buncher[1]

    #delimit ;
    twoway (scatter bunch_rel ebitda2atsd_bin if early == 1) 
        (lfit bunch_rel ebitda2atsd_bin if early == 1, lc(red) lp(dash))
        (scatter bunch_rel ebitda2atsd_bin if early == 0, ms(S) mc(maroon))
        (lfit bunch_rel ebitda2atsd_bin if early == 0, lc(red) lp(dash)),
        ytitle("P(Bunch at Threshold, Indexed in % to Bin 1)")
        yscale(range(0 250)) ylabel(0 (50) 250)
        xtitle("Bin of EBITDA Volatility") `graphopts'
        legend(order(1 "Early" 3 "Late"));
    graph export "ebitda2atsd_bunch.pdf", as(pdf)  replace;

    gen ebitda2atsd_bin2 = ebitda2atsd_bin - 0.5 if early == 1;
    replace ebitda2atsd_bin2 = ebitda2atsd_bin + 0.5 if early == 0;

    twoway
        (bar variance_negshare ebitda2atsd_bin2 if early == 1, color(navy) barwidth(1))
        (bar variance_negshare ebitda2atsd_bin2 if early == 0, color(maroon) barwidth(1)),
        ytitle("Variance from Loss (%)")
        yscale(range(0 100)) ylabel(0(50)100)
        xtitle("Bin of EBITDA Volatility") $gpr legend(off)
        ysize(2) scale(2);
    graph export "ebitda2atsd_loss_bunch.pdf", replace;
    #delimit cr

    ****************************************************************************
    * 2. Autocorrelation graph
    ****************************************************************************
    tempfile early
    load_analysis_data bunch_early
    gen early = 1
    save `early'

    load_analysis_data bunch_late
    gen early = 0
    append using `early'

    #delimit
    twoway 
        (rcap beta_u beta_l lag if early == 1)
        (connect beta lag if early == 1, msym(o) mc(navy) lc(navy) lp(dash))
        (rcap beta_u beta_l lag if early == 0, mc(maroon) lc(maroon))
        (connect beta lag if early == 0, msym(s) mc(maroon) lc(maroon) lp(dash))
        , 
        yline(0, lc(red)) 
        xsca(range(-.5 7.5)) xlab(0(1)7)
        ysca(range(-.12 .01)) ylab(-.12(.02)0)
        ytitle("P(Future Bunch | Bunch in Year 0)") 
        xtitle("Years from Bunch") $gpr $graphconfig 
        legend(order(2 "Early" 4 "Late"));
    #delimit cr
    graph export "prob_future_bunch.pdf", replace

    ****************************************************************************
    * 3. Cumulative graph
    ****************************************************************************
    
    /* Output from logfile generated with internal code file.
    * Hard-coded for the graph here to make formatting easier.
    * bunchstack == 1
    -------------+----------------------------------------------------------------
    zero |          0  (omitted)
    d0 |   .0134267   .0039007     3.44   0.001     .0057773    .0210761
    f1 |   .0297056   .0041395     7.18   0.000      .021588    .0378232
    f2 |   .0354237   .0046481     7.62   0.000     .0263087    .0445387
    f3 |   .0425676   .0047986     8.87   0.000     .0331575    .0519776
    early |   .0127871   .0168759     0.76   0.449    -.0203067    .0458809
    ebitda2at |   .0002416   .0002021     1.20   0.232    -.0001547     .000638
    _cons |   .0709047   .0124548     5.69   0.000     .0464807    .0953287
    ------------------------------------------------------------------------------
    * bunchstack == 0 => Within 2.5K for early and 25K for late
    -------------+----------------------------------------------------------------
    zero |          0  (omitted)
    d0 |   .0187498   .0012897    14.54   0.000     .0162219    .0212778
    f1 |   .0212709   .0013381    15.90   0.000     .0186481    .0238936
    f2 |   .0226152   .0014148    15.99   0.000     .0198421    .0253883
    f3 |   .0196752   .0014555    13.52   0.000     .0168223     .022528
    early |   .0229998   .0041409     5.55   0.000     .0148833    .0311163
    ebitda2at |  -2.75e-08   9.97e-08    -0.28   0.782    -2.23e-07    1.68e-07
    _cons |   .0864324    .001756    49.22   0.000     .0829905    .0898744
    ------------------------------------------------------------------------------
    */

    clear
    set obs 5
    gen event = _n - 2
    gen coef_bunch = 0 if event == -1
    replace coef_bunch = 0.0134267 if event == 0
    replace coef_bunch = 0.0297056 if event == 1
    replace coef_bunch = 0.0354237 if event == 2
    replace coef_bunch = 0.0425676 if event == 3

    gen se_bunch = 0 if event == -1
    replace se_bunch = 0.0039007 if event == 0
    replace se_bunch = 0.0041395 if event == 1
    replace se_bunch = 0.0046481 if event == 2
    replace se_bunch = 0.0047986 if event == 3

    gen coef_bunch_upper = coef_bunch + 1.96 * se_bunch
    gen coef_bunch_lower = coef_bunch - 1.96 * se_bunch

    gen coef_nearbunch = 0 if event == -1
    replace coef_nearbunch = 0.0187498 if event == 0
    replace coef_nearbunch = 0.0212709 if event == 1
    replace coef_nearbunch = 0.0226152 if event == 2
    replace coef_nearbunch = 0.0196752 if event == 3

    gen se_nearbunch = 0 if event == -1
    replace se_nearbunch = 0.0012897 if event == 0
    replace se_nearbunch = 0.0013381 if event == 1
    replace se_nearbunch = 0.0014148 if event == 2
    replace se_nearbunch = 0.0014555 if event == 3

    gen coef_nearbunch_upper = coef_nearbunch + 1.96 * se_nearbunch
    gen coef_nearbunch_lower = coef_nearbunch - 1.96 * se_nearbunch

    #delimit ;
    twoway 
        (connect coef_bunch event, lp(dash))
        (rcapsym coef_bunch_upper coef_bunch_lower event, color(navy) msym(none))
        (connect coef_nearbunch event, color(maroon))
        (rcapsym coef_nearbunch_upper coef_nearbunch_lower event, color(maroon) msym(none)),
        yline(0) 
        xlab(-1 "t-1" 0 "t0" 1 "t+1" 2 "t+2" 3 "t+3")
        ylab(, format(%9.2f))
        ytitle("Investment/Initial Capital")
        xtitle("Event Year")
        graphregion(color(white) lwidth(large)) bgcolor(white) 
        legend(order(1 "Bunchers" 3 "Near Bunchers"));
    graph export "cumulative_investment_bunch.pdf", replace;
    #delimit cr

end

 /* %> */



****************************************************
* Figure 6 panel A: Fiscal Q4 Spikes in Model Simulated Data

capture program drop fig_model_agg_spikes /* %< */
program define fig_model_agg_spikes

	load_analysis_data model

	tempfile orig
	save `orig'
	
	foreach version in ov_beta975 hycta_beta975 full_beta975 {
	    
		use `orig', clear
		keep if version == "`version'"

		* Only keep first 25 periods and investment data
		drop if period > 25
		drop gos ktax k 
		capture drop kbar g

		* Calculate average invratio for each quarter
		collapse (mean) inv, by(year quarter period)

		* Calculate average yearly investment for each year
		egen avg_year_investment = mean(inv), by(year)
	
		* Calculate aggregate invratio for each period
		gen invratio = 100 * (inv / avg_year_investment)

		* Rename variables to make it clear which version they're from
		ren invratio invratio_`version'
		
		* Create a point showing the investment in quarter 4 for each year
		gen invratio_q4_`version' = invratio_`version' if quarter == 4

		* Make variable names lowercase
		rename *, lower
				
		tempfile agg_spike_`version'
		save "`agg_spike_`version''"
	}
	
	* Merge the datasets together 
	use "`agg_spike_ov_beta975'", clear
	merge 1:1 year quarter using "`agg_spike_hycta_beta975'", assert(3) nogen
	merge 1:1 year quarter using "`agg_spike_full_beta975'", assert(3) nogen
	
	* Make graph
	#delimit ;
	graph twoway
		(line invratio_full_beta975 period, lcolor(dknavy) lpattern(solid) lwidth(medthick))
		(scatter invratio_q4_full_beta975 period, mcolor(dknavy) msize(small) msymbol(s))
		(line invratio_hycta_beta975 period, lcolor(forest_green) lpattern(dash))
		(scatter invratio_q4_hycta_beta975 period, mcolor(forest_green) msize(small) msymbol(o))
		(line invratio_ov_beta975 period, lcolor(cranberry) lpattern(dash_dot))
		(scatter invratio_q4_ov_beta975 period, mcolor(cranberry) msize(small) msymbol(t)),
		graphregion(color(white) lwidth(large)) bgcolor(white)
		legend(order(1 "Full Model" 3 "Depreciation Motive" 5 "Baseline") rows(1))
		ytitle("Aggregate Investment Ratio") xtitle("Period");
	#delimit cr
	graph export "agg_spike_graph.pdf", replace
	
end /* %> */

*************************************************************************
* Figure 6 Panels B,C,D : Depreciation versus Option Value Motives in Model Simulations

capture program drop fig_model_cumlcompare /*%<*/
program define fig_model_cumlcompare

	model_bootstrap
	
    use "$intermediate/simulation_bootstrap_20200517.dta", clear

    ttest spike_coeff_full = spike_coeff_hycta if event_time == 12

    * Pick a simulation where the diff in t = 12 coefficients is close to the
    * mean reported in the text and the initial spikes are similar. 
    * Goal of the figure is to show a representative simulation. After applying
    * the above restriction, we just grab one of the simulations by index
    * number, which is hard-coded for replicability.
    * Moments of the distribution are reported in the text.

    gen diff = spike_coeff_full - spike_coeff_hycta
    sum diff
    gen near_mean = abs(diff - `r(mean)') < .01 & diff > 0
    gen seriesholder = near_mean * (event_time == 12)
    bys iteration: egen keeper = max(seriesholder)
    keep if keeper == 1
    keep if iteration == 456 
    drop if event_time <= 0
    sort event_time

    replace spike_coeff_full = 100 * spike_coeff_full
    replace spike_coeff_hycta = 100 * spike_coeff_hycta
    replace spike_coeff_ov = 100 * spike_coeff_ov

    #delimit ;
    twoway 
        (connect spike_coeff_full event_time, msym(s) mc(navy) lc(navy) lp(solid))
        (connect spike_coeff_hycta event_time, msym(o) mc(forest_green) 
			lc(forest_green) lp(dash))
        (connect spike_coeff_ov event_time, msym(t) mc(cranberry) lc(cranberry)
			lp(dash_dot)),
        yline(0, lc(gray)) 
        xsca(range(0.5 12.5)) xlab(1(1)12)
        ytitle("Cumulative Avg. Investment Relative to Baseline (%)") 
        xtitle("Quarter") $graphconfig
        xline(4, lp(dash) lc(gray))
        legend(order(1 "Full Model" 2 "Depreciation Motive" 3 "Baseline") col(3));
    #delimit cr
    graph export "cum_model_comparison_vMean.pdf", replace

end/*%>*/

capture program drop fig_model_spikes_by_earn_bar /* %< */
program define fig_model_spikes_by_earn_bar

	* Figures 6c and 6d
	
	syntax, bins(integer) 
	
	* Only keep a subset of firms - all firms are identical and this makes code faster 
	keep if firm_no <= 500
	
	tempfile orig
	save `orig'

	clear
	set obs `bins'
	gen n = _n
	
	tempfile spikes_by_earn_mean spikes_by_earn_vol
	save "`spikes_by_earn_mean'", replace 
	save "`spikes_by_earn_vol'", replace 

	foreach version in ov_beta975 hycta_beta975 full_beta975 {
		
		use `orig', clear
		keep if version == "`version'"
		
		*** Create variables needed for graph ***
		
		* Invratio
		forval i = 1/4 {
			egen q`i'_invest = max(inv) if quarter == `i', by(firm_no year)
			egen q`i'_investment = max(q`i'_invest), by(firm_no year)
			drop q`i'_invest
		}
		
		egen year_investment = mean(inv), by(firm_no year)
		gen invratio = q4_investment/ year_investment
		
		* Annual gross operating surplus is the sum of g and gos in quarter 4
		* in full and hyc model
		if strpos("`version'", "hyc") > 0 | strpos("`version'", "full") > 0 {
			gen tot_gos1 = gos + g if quarter == 4
			egen tot_gos = max(tot_gos1), by(firm_no year)
			drop tot_gos1 q4_invest 
		}
		else {
			* In the baseline model annual gos is just sum of quarterly gos q1-q4
		    egen tot_gos = sum(gos), by(firm_no year)
		}
				
		* Collapse over quarters
		collapse (last) k ktax (max) invratio tot_gos, by(firm_no year)
			
		* Censor invratio similar to in empirical analysis 
		* and then multiply invratio by 100 to get a percentage
		
		drop if invratio > 5
		replace invratio = invratio * 100
				
		* Calculate mean and st dev of ebitda/assets
		egen mean_cap = mean(k), by(firm_no)
		gen earn_mean = tot_gos / mean_cap
		gen gos_cap_ratio = tot_gos / k
		egen earn_vol = sd(gos_cap_ratio), by(firm_no)
	
		gen n = _n
		
		tempfile pre_collapse_`version'
		save `pre_collapse_`version''		
	}
	
	* Calculate bins of within-firm mean and st dev of ebitda/assets for the full model
	use `pre_collapse_full_beta975', clear
		
	preserve 
		foreach var of varlist earn_vol earn_mean {
			pctile `var'_cut = `var', nq(10)
		}
		
		keep *_cut n
		
		tempfile cuts
		save `cuts'
		
	restore
	
	foreach version in ov_beta975 hycta_beta975 full_beta975 {
		
		use `pre_collapse_`version'', clear
		merge 1:1 n using `cuts', keep(1 3)
		
		drop n _merge
		
		* Create quantiles with cuts from the full model
		foreach var of varlist earn_vol earn_mean {
			xtile `var'_bin = `var', cut(`var'_cut)
		}
		* Rename to make it obvious which version the variables are from
		ren (invratio earn_mean_bin earn_vol_bin) (invratio_`version' ///
			mean_bin_`version' vol_bin_`version')
		ren (earn_mean earn_vol) (mean_`version' vol_`version')
		
		* Collapse to get average invratio - 100 by bin
		foreach type in mean vol {
			preserve
				collapse (mean) invratio_`version' `type'_`version', by(`type'_bin_`version')
				replace invratio_`version' = invratio_`version' - 100
				tempfile `data'_`type'
				save ``data'_`type'', replace
			restore
		}
		
	* Merge together
		foreach type in mean vol {
			use ``data'_`type'', clear
			gen n = _n
			merge 1:1 n using "`spikes_by_earn_`type''", nogen 
			save "`spikes_by_earn_`type''", replace
		}
	}

	* Make graphs	
	use "`spikes_by_earn_mean'", clear
	
	tostring mean_full_beta975, format("%3.2f") gen(earn_mean) force

	* Figure 6c	
	
	#delimit ;
		graph bar (mean) invratio_full_beta975 invratio_hycta_beta975 invratio_ov_beta975, 
		over(earn_mean, sort(n)) 
		legend(label(1 "Full Model") label(2 "Depreciation Motive") 
			   label(3 "Baseline") rows(1))
		ytitle("Mean Investment Ratio") b1title("GOS / Mean Assets") 
		yscale(r(-20 40)) 
		ylabel(-20 "80" -10 "90" 0 "100" 10 "110" 20 "120" 30 "130" 40 "140")
		bar(1, color(dknavy)) bar(2, color(forest_green)) bar(3, color(cranberry))
		$graphconfig;
	#delimit cr 
	graph export "spikes_by_earn_mean_`bins'bins_trim5.pdf", replace

	* Figure 6d 
	
	use "`spikes_by_earn_vol'", clear
	tostring vol_full_beta975, format(%3.2f) gen(earn_vol) force
	#delimit ;
		graph bar (mean) invratio_full_beta975 invratio_hycta_beta975 invratio_ov_beta975, 
		over(earn_vol, sort(n)) 
		legend(label(1 "Full Model") label(2 "Depreciation Motive") 
		       label(3 "Baseline") rows(1))
		ytitle("Mean Investment Ratio") b1title("SD(GOS / Assets)")
		yscale(r(0 15)) ylabel(0 "100" 5 "105" 10 "110" 15 "115")
		bar(1, color(dknavy)) bar(2, color(forest_green)) bar(3, color(cranberry))
		$graphconfig;
	#delimit cr 
	graph export "spikes_by_earn_vol_`bins'bins_trim5.pdf", replace
	
end /* %> */

**********************************************************************
* Figure 7: Spikes in Capital Goods Shipments, Inventories, and Prices

capture program drop fig_seasonal_agg /* %< */
program define fig_seasonal_agg

    collapse  vs_n ti_n uo_n no_n , by(industry month)

	* Figure 7a: Spikes in Capital Goods Shipments (1958–2016)
	
    #delimit;
    twoway 
		(connect vs_n month if industry=="NDE", lp(dash) ms(O) msize(medlarge)) 
		(connect vs_n month if industry=="COG", lp(dash) ms(S) msize(medlarge)), 
		xsc(range(1 12)) xlab(1(1)12) yline(100, lc(gray)) ysc(range(80 120))  
		ylab(80 (10) 120) ysize(3) 
		graphregion(color(white) lwidth(large)) bgcolor(white)
		ytitle("Monthly Shipments/Average Shipments (%)")
		xtitle("Calendar Month") 
		legend(order (1 "Capital Goods Shipment Value"
				      2 "Consumer Goods Shipment Value"));
    #delimit cr
    graph export "M3combine.pdf", as(pdf) replace
    
	* Figure 7b: Spikes in Capital Goods Shipments and Inventories (1958–2016)

	#delimit;
	twoway 
		(connect vs_n month if industry=="NDE", ms(O) msize(medlarge) lp(dash)) 
		(connect ti_n month if industry=="NDE", ms(S) msize(medlarge) lp(dash) yaxis(2)),
		xtitle(" Calendar Month") ytitle("Monthly  Shipment Value to Average Value (%)")   
	    ytitle("Monthly  Inventory to Average  Inventory (%)", axis(2))  
		graphregion(color(white) lwidth(large)) bgcolor(white) 
		xlabel(1(1) 12)  ysize(3) yscale(r(80 120)) ysc(r(80 120) axis(2)) 
		ylabel(80(10)120) ylabel(80(10)120, axis(2)) yline(100, lc(black))
		legend(order (1 "Capital Goods Shipment Value" 
					  2 "Capital Goods Total Inventory" ));	
	#delimit cr
	graph export "M3inventory.pdf", as(pdf)  replace 

end /* %> */

capture program drop fig_seasonal_agg_prices /* %< */
program define fig_seasonal_agg_prices

	* Figure 7c: Spikes in Capital Goods Shipments and Prices (1998–2016)

	#delimit;
	twoway	
		(connect vs_ratio_NDE month, ms(O) msize(medlarge) lp(dash)) 
		(connect ppi_ratio_vw month, ms(S)  msize(medlarge) lp(dash) yaxis(2) ), 
		xtitle(" Calendar Month") 
		ytitle("Monthly  Shipment Value to Average Value (%)")  
		ytitle("Monthly PPI to Averge PPI (%)", axis(2))  
		legend(order (1 "Capital Goods Shipment Value" 
					  2 "Capital Goods Producer Price Index" ))   
		graphregion(color(white) lwidth(large)) bgcolor(white) 
		xlabel(1(1) 12)  ysize(3) yscale(r(80 120)) yscale(r(80 120) axis(2)) 
		ylabel(80(10)120, axis(2)) ylabel(80(10)120) yline(100, lc(black));
	#delimit cr		
	graph export "M3vs_ppi.pdf", as(pdf)  replace 
	
end /* %> */

*************************************************************
* Figure 8: Supplier Q4 Inventory and Sales Spikes 

capture program drop fig_supplier /* %< */
program define fig_supplier

	* ff12: 3 manufacturing 6 business equipment
	keep if ff12 == 3 | ff12 == 6
	order gvkey fyear sale c_gvkey customer_sale invtq* c_capx43 c_capx*
	drop if invtq1 == . | invtq2 == . | invtq3 == . | invtq4 == .
	drop if invtq1 == 0 & invtq2 == 0 & invtq3 == 0 & invtq4 == 0
	drop if c_capx43 == .
	winsor2 c_capx43, replace cut(3 97)
	replace c_capx43 = 5 if c_capx43 > 5
	su c_capx43, de
	gen capx43_ratio = c_capx43 * customer_sale / sale
	bysort gvkey fyear: egen total_csale = sum(customer_sale)
	gen ratio = total_csale/sale
	bysort gvkey: egen minratio = min(ratio)
	keep if minratio > 0.1
	
	* Trace out suppliers' inventory seasonality within a year
	gen invbase14 = (invtq1 + invtq2 + invtq3 + invtq4)/4
	gen invbase13 = (invtq1 + invtq2 + invtq3)/3
	forvalues i=1/4 {
		gen inv`i'_4 = invtq`i'/invbase14
		winsor2 inv`i'_4, replace cut(1 99)
		su inv`i'_4, de
	} 
	gen inv4_3 = invtq4/invbase13
	winsor2 inv4_3, replace cut(1 99)
	su inv4_3, de

	collapse inv* c_capx43 sale4_3 sale4_4 (sum) capx43_ratio, by(gvkey fyear)
	 
	foreach x in inv4_3 inv3_4 inv4_4 sale4_4 sale4_3 c_capx43 capx43_ratio {
		winsor2 `x', replace cut(1 99)
		replace `x' = `x' * 100
	}
		 
	* Figure 8a
	
	binscatter sale4_3 c_capx43, absorb(gvkey) nq(40) name(saleew, replace) ///
	ytitle("Supplier Sales Q4/Ave(Q1-Q3) %") ///
	xtitle("Customer CAPEX Q4/Ave(Q1-Q3) % ") ///
	graphregion(color(white) lwidth(large)) bgcolor(white) ///
	legend(off)  ysize(3.5)  yscale(r(100 120)) ylabel(100(5)120)
	graph export "suppliersale_ew.pdf", as(pdf)  replace
	
	* Figure 8b
	
	binscatter inv4_3 c_capx43, absorb(gvkey) nq(40) name(invew, replace) ///
	ytitle( "Supplier Inventory Q4/Ave(Q1-Q3) %") ///
	xtitle("Customer CAPEX Q4/Ave(Q1-Q3) % ") ///
	graphregion(color(white) lwidth(large)) bgcolor(white) ///
	legend(off)  ysize(3.5)
	graph export "supplierinv_ew.pdf", as(pdf)  replace
	
end /* %> */

*****************************************************************
* Figure 9: Capital Lending Volume and Interest Rate (2005–2018) 

capture program drop fig_lending /* %< */
program define fig_lending

	* Figure 9a
	
    gen time = ym(year,month)
    gen time2 = time
    format time %tm 

    forv j=551(12)707 {
        local spikes `spikes', `j'
    }

    #delimit ;
    twoway 
		(line total_volumebillionsofusd time) (scatter total_volumebillionsofusd
		time if inlist(time `spikes'), mcolor(red)),
        ytitle("New Business Volume (Loans & Leases, $ Billions)") xtitle("Month")
        xlabel(551 (24) 707) legend(off) xsize(7)
        graphregion(color(white) lwidth(large)) bgcolor(white);
    #delimit cr
    graph export "newbusiness.pdf", as(pdf) replace

end /* %> */

capture program drop fig_loanrate /* %< */
program define fig_loanrate
	
	* Figure 9b
	
	#delimit ;
	twoway	
		(connect commequp250k month, ms(O) msize(medium) lp(solid)) 
		(connect commreal1mil month, ms(D)  msize(medium) lp(dash))
		(connect personalunsec month, ms(S)  msize(medium) lp(shortdash)),
		xtitle(" Calendar Month") ytitle("Monthly Rate to Average  Rate (%)") 
		legend(order (1 "Commercial Equipment" 
					  2 "Commercial Real Estate" 
					  3 "Personal Loan")  cols(3)) 
		graphregion(color(white) lwidth(large)) bgcolor(white) 
		xlabel(1(1) 12) ysize(3) yscale(r(90 110)) ylabel(90(2)110) 
		yline(100, lc(black));
	#delimit cr		
	graph export "loanratemean.pdf", as(pdf)  replace 

end /* %> */

*****************************
* Table 2: Summary Statistics 

capture program drop table_sumstats /* %< */
program define table_sumstats

	* Table 2a
	
    tempfile orig
    save `orig'
	replace capx4_3 = capx4_3/100

    drop if at == . | at == 0
    la var cf2at "Cash Flow/Assets"
    la var res2at "Cash/Assets"
    la var ebitda2at "EBITDA/Assets"
    la var capx4_3 "CAPEX 4/3 (\%)"
    la var sale4_3 "Sales 4/3 (\%)"
    la var capx2ppe "CAPEX/PPE"
    la var ppe "PPE (\\$ M)"
    la var sale "Sales (\\$ M)"
    la var at "Assets (\\$ M)"
    la var capx "CAPEX (\\$ M)"
    la var dp "Depreciation (\\$ M)"
    la var mkt2bk "M/B"

    #delimit ;
    estpost tabstat at dp capx ppe sale  mkt2bk cf2at res2at ebitda2at
        capx2ppe capx4_3 sale4_3, statistics(n mean median sd p10 p90)
        columns(stats);

    esttab . , cell("count(fmt(%9.0fc)label(N)) mean(fmt(%9.2fc)label(Mean))
        p50(fmt(%9.2fc)label(Median)) sd(fmt(%9.2fc)label(SD))
        p10(fmt(%9.2fc)label(P10)) p90(fmt(%9.2fc)label(P90)) ")
    varlabels(`e(labels)')  label noobs var(23) nonumb ;

    esttab . using "summary.tex", booktabs  nomti
        cell(" count(fmt(%9.0fc)label(N)) mean(fmt(%9.2fc)label(Mean))
        p50(fmt(%9.2fc)label(Median)) sd(fmt(%9.2fc)label(SD))
        p10(fmt(%9.2fc)label(P10)) p90(fmt(%9.2fc)label(P90))")
    varlabels(`e(labels)') label noobs var(23) nonumb replace;
    #delimit cr

    use `orig', clear

end /* %> */

capture program drop table_sumstats_intl /* %< */
program define table_sumstats_intl

	* Table 2b
	
    drop if at == . | at == 0
    la var cf2at "Cash Flow/Assets"
    la var res2at "Cash/Assets"
    la var ebitda2at "EBITDA/Assets"
    la var capx4_3 "CAPEX 4/3 (\%)"
    la var sale4_3 "Sales 4/3 (\%)"
    la var capx2ppe "CAPEX/PPE"
    la var ppent "Property, Plant and Equipment (Mils)"
    la var sale "Sales (Mils)"
    la var at "Assets (Mils)"
    la var capxy4 "Capital Expenditures (Mils)"
    la var dpy4 "Depreciation and Amortization (Mils)"
    la var mkt2bk "M/B"

    estpost tabstat mkt2bk cf2at res2at ebitda2at capx2ppe capx4_3 sale4_3, ///
		    statistics(n mean median sd p10 p90) columns(stats)

    esttab . , cell("count(fmt(%9.0fc)label(N)) mean(fmt(%9.2fc)label(Mean)) p50(fmt(%9.2fc)label(Median)) sd(fmt(%9.2fc)label(SD)) p10(fmt(%9.2fc)label(P10)) p90(fmt(%9.2fc)label(P90)) ") ///
	    varlabels(`e(labels)')  label noobs var(23) nonumb 

    esttab . using "internationalsummary.tex", booktabs  nomti ///
		cell(" count(fmt(%9.0fc)label(N)) mean(fmt(%9.2fc)label(Mean)) p50(fmt(%9.2fc)label(Median)) sd(fmt(%9.2fc)label(SD)) p10(fmt(%9.2fc)label(P10)) p90(fmt(%9.2fc)label(P90))") ///
	    varlabels(`e(labels)') label noobs var(23) nonumb replace

end /* %> */

* Table 3: Fiscal Q4 CAPEX Spikes and Tax Status 
* Note: This table uses internal data 

******************************************************
* Table 4: Investment Spikes and Financial Constraints

capture program drop table_constraints /* %< */
program define table_constraints

    tempfile orig
    save `orig'
	replace capx4_3 = capx4_3/100

    drop if fyear == 1983
    gen d8487 = (inlist(fyear, 1984, 1985, 1986, 1987))
    la var d8487 "D(84-87)"
    gen nodiv = 1 if paydiv == 0
    replace nodiv = 0 if paydiv == 1
	
    foreach x of varlist nodiv junkrating fp fp2 {
        gen inter`x' = d8487 * `x'
        la var inter`x' "D(1984-1987)*`x'" 
        eststo tra`x': qui reghdfe capx4_3 d8487 `x' inter`x' $controls1, ///
			   absorb(gvkey) cluster(gvkey)
        qui: estadd local firm "Yes", replace
        qui: estadd local controls "Yes", replace
        qui: estadd local year "No", replace
    }
	gen interl_at = d8487 * l_at
	la var interl_at "D(1984-1987)*ln(assets)" 
	eststo tral_at: qui reghdfe capx4_3 d8487 interl_at $controls1, ///
		   absorb(gvkey) cluster(gvkey)
	qui: estadd local firm "Yes", replace
	qui: estadd local controls "Yes", replace
	qui: estadd local year "No", replace
    
	esttab  tral_at tranodiv trajunkrating trafp trafp2 using "fc.tex", ///
	scalar("N Observations" "r2_a Adj R-Squared" "controls Controls" ///
	       "year Year FE" "firm Firm FE" ) ///
	keep(d8487 inter*) order(d8487 inter*)  ar2(%8.2f) ///
	cells(b(star fmt(%9.2f)) se(par) ) starlevels(* .10 ** .05 *** .01 )  ///
	style(tex) replace  ///
	label collabels(none) varlabels(_cons Constant)  nomti nodep
	
    use `orig', clear

end /* %> */

***********************************************
* Table 5: Cumulative Effect of Q4 CAPEX Spikes
* Note: This table uses internal data

**********************************************************************
* Table 6: Investment Spikes and Complicated Firms: Use it or Lose it?

capture program drop table_budget /* %< */
program define table_budget

use $analysis/capxannual_trim5.dta, clear
merge 1:1 gvkey fyear using $analysis/orbis_layer.dta, keep(1 3) nogen
merge 1:1 gvkey fyear using $analysis/Compustat/exeownership.dta, keep(1 3) nogen
la var maxlayer "\# Layers"
winsor2 maxlayer top5exe_pct ceo_pct cfo_pct allexe_pct, cut(1 99) replace 

foreach x of varlist maxlayer  top5exe_pct ceo_pct cfo_pct allexe_pct {
egen `x'sd=sd(`x')
gen `x'_sd=`x'/`x'sd
replace `x'=`x'_sd
}

 foreach x of varlist  maxlayer  allexe_pct { 
eststo `x': qui reghdfe capx4_3 `x' $controls1 , absorb(fyear) cluster(gvkey)
qui: estadd local year "Yes"
qui: estadd local controls "Yes"
 }
 
	order gvkey capx total_capxs sale total_sales at total_ias 
		
	* Require that the total capx and sales cover at least 80% of the numbers
	* reported on Compustat
	drop if total_capxs / capx < 0.8
	drop if total_sales / sale < 0.8
	gen seg_sort = (seg_num > 1 & seg_num != .)
	gen sic3_sort = (sic3_count > 1 & seg_num != .)
	gen sic2_sort = (sic2_count > 1 & seg_num != .)
	la var seg_sales_sd2mean " $\frac{SD(Sales)}{Ave(Sales)}$"
	la var seg_capxs_sd2mean " $\frac{SD(CAPEX)}{Ave(CAPEX)}$ "
	la var seg_num "\# Segments"
	la var sic2_count "\# SIC2"
	la var allexe_pct "Exec Own \%"

	foreach x of varlist seg_num sic3_count sic2_count {
		egen `x'sd = sd(`x')
		gen `x'_sd = `x' / `x'sd
		replace `x' = `x'_sd
	}

	foreach x of varlist seg_num sic2_count {
		eststo `x': qui reghdfe capx4_3 `x' $controls1 , ///
						absorb(fyear) cluster(gvkey)
		qui: estadd local year "Yes"
		qui: estadd local controls "Yes"
		eststo `x'4: qui reghdfe capx4_4 `x' $controls1 , ///
						 absorb(fyear) cluster(gvkey)
		qui: estadd local year "Yes"
		qui: estadd local controls "Yes"
		eststo `x'4y: qui reghdfe capx4_y `x' $controls1 , ///
						  absorb(fyear) cluster(gvkey)
		qui: estadd local year "Yes"
		qui: estadd local controls "Yes"
	}
	 
	esttab  seg_num sic2_count maxlayer allexe_pct using "budget.tex",   ///
	scalar("N Observations" "r2_a Adj R-Squared" "controls Controls"  ///
		   "year Year FE" )   ///
	keep( seg_num sic2_count  maxlayer allexe_pct  ) ///
	cells(b(star fmt(%9.1f)) se(par) )   ar2(%8.2f) ///
	starlevels(* .10 ** .05 *** .01 ) style(tex) replace  ///
	label collabels(none) varlabels(_cons Constant) nodep nomti	

end /* %> */

**********************************************************************
* Table 8: The Effect of a Corporate Tax Cut on Investment in Model Simulations

capture program drop table_model_elasticity/*%<*/
program define table_model_elasticity

	/* This program calculates the results on Table VII. Performs a DiD using a tax-policy experiment.
	10.000 firms for 150 years. At year = 75 there is a tax change on half of firms from 35% to 21%
	*/

    local deltahat = 0.119

    foreach k of varlist k* {
        replace `k' = `k'[_n-1] if firm==firm[_n-1] & year==year[_n-1] & quarter==2
        replace `k' = `k'[_n-2] if firm==firm[_n-2] & year==year[_n-2] & quarter==3
        replace `k' = `k'[_n-3] if firm==firm[_n-3] & year==year[_n-3] & quarter==4
    }

    collapse (sum) inv* (mean) k* (sum) gos*, by(firm year group post)

	by firm: egen kpre = total(cond(year == 74, k, .))   

    xtset firm year

    gen inv_rate = inv / kpre
    gen loss = gos - 4 * `deltahat' * kbar < 0

    * calibrate z = 0.881 from Zwick Mahon
    gen taxterm_21 	= (1-0.21*0.881)/(1-0.21)
    gen taxterm_35 	= (1-0.35*0.881)/(1-0.35)

    gen taxterm = taxterm_35
    replace taxterm = taxterm_21 if group == 1 & post == 1

    gen log_inv = log(inv)

    gen keeprate_21 = log(1-0.21)
    gen keeprate_35 = log(1-0.35)
    
    gen keeprate = keeprate_35
    replace keeprate = keeprate_21 if group == 1 & post == 1

    keep if inrange(year, 74, 77)

    * 1. Diff-in-diff table
    preserve
    collapse (mean) inv inv_rate keeprate taxterm, by(group post loss)
    sort group loss post 
    list group loss post inv_rate

    sort group loss post 
    by group loss: gen inv_rate_pre = inv_rate[1]
    by group loss: gen inv_rate_post = inv_rate[2]
    gen inv_rate_diff = inv_rate_post - inv_rate_pre

    sort loss post group
    by loss post: gen time_diff = inv_rate[2] - inv_rate[1]
    sort loss group post
    by loss group: gen time_diff_pre = time_diff[1]
    by loss group: gen time_diff_post = time_diff[2]
    gen diff_in_diff = time_diff_post - time_diff_pre

    gen semi_elasticity = 100 * diff_in_diff / inv_rate_pre
    sum semi_elasticity if group == 1 & loss == 0 & post == 1
    replace semi_elasticity = round(semi_elasticity, 0.001)
    sum semi_elasticity if group == 1 & loss == 0 & post == 1

    gen taxterm_diff = (1 - 0.21 * 0.881)/(1-0.21) - (1 - 0.35 * 0.881)/(1 - 0.35) 
    gen taxterm_elasticity = diff_in_diff / taxterm_diff
    replace taxterm_elasticity = round(taxterm_elasticity, 0.001)

    format *_elasticity %04.3f

    foreach v of varlist inv_rate_* time_diff* diff_in_diff {
        replace `v' = round(100 * `v', 0.01)
        format `v' %4.2f
    }

    gen diff_name = "Diff"
    gen group_name = cond(group == 1, "Tax Rate Cut", "No Tax Change")

    gen tab = "\begin{tabular}{lccc}"
    gen toprule = "\toprule"
    gen midrule = "\midrule"
    gen cmidrule = "\cmidrule(lr){2-4}"
    gen botrule = "\bottomrule"
    gen end = "\end{tabular}" 

    gen c = " "

    gen title1 = "& I/K Pre (\%) & I/K Post (\%) & Diff"
    gen linespace = "\addlinespace"
    gen semi_name = "Semi-elasticity (\% of I/K Pre)"
    gen elas_name = "Tax Term Coefficient"

    local tabname = "model_did_nonloss20240101"
    * Output
    listtex tab if _n == 1 using "`tabname'.tex", replace rstyle(none)
    listtex toprule if _n == 1, appendto("`tabname'.tex") rstyle(none)
    listtex title1 if _n == 1, appendto("`tabname'.tex") rstyle(tabular)
    listtex cmidrule if _n==1 , appendto("`tabname'.tex") rstyle(none)    

    listtex group_name inv_rate_pre inv_rate_post inv_rate_diff if group == 0 & loss == 0 & post == 1, ///
            appendto("`tabname'.tex") rstyle(tabular)
    listtex group_name inv_rate_pre inv_rate_post inv_rate_diff if group == 1 & loss == 0 & post == 1, ///
            appendto("`tabname'.tex") rstyle(tabular)
    listtex midrule if _n==1 , appendto("`tabname'.tex") rstyle(none)    
    listtex diff_name time_diff_pre time_diff_post diff_in_diff if group == 0 & loss == 0 & post == 1, ///
            appendto("`tabname'.tex") rstyle(tabular)
    listtex linespace if _n == 1, appendto("`tabname'.tex") rstyle(none)
    listtex semi_name c c semi_elasticity if group == 1 & loss == 0 & post == 1, appendto("`tabname'.tex") rstyle(tabular)
    listtex elas_name c c taxterm_elasticity if group == 1 & loss == 0 & post == 1, appendto("`tabname'.tex") rstyle(tabular)
    listtex botrule if _n == 1, appendto("`tabname'.tex") rstyle(none)  
    listtex end if _n == 1, appendto("`tabname'.tex") rstyle(none)

    local tabname = "model_did_loss20240101"
    * Output
    listtex tab if _n == 1 using "`tabname'.tex", replace rstyle(none)
    listtex toprule if _n == 1, appendto("`tabname'.tex") rstyle(none)
    listtex title1 if _n == 1, appendto("`tabname'.tex") rstyle(tabular)
    listtex cmidrule if _n==1 , appendto("`tabname'.tex") rstyle(none)    

    listtex group_name inv_rate_pre inv_rate_post inv_rate_diff if group == 0 & loss == 1 & post == 1, ///
            appendto("`tabname'.tex") rstyle(tabular)
    listtex group_name inv_rate_pre inv_rate_post inv_rate_diff if group == 1 & loss == 1 & post == 1, ///
            appendto("`tabname'.tex") rstyle(tabular)
    listtex midrule if _n==1 , appendto("`tabname'.tex") rstyle(none)    
    listtex diff_name time_diff_pre time_diff_post diff_in_diff if group == 0 & loss == 1 & post == 1, ///
            appendto("`tabname'.tex") rstyle(tabular)
    listtex linespace if _n == 1, appendto("`tabname'.tex") rstyle(none)
    listtex semi_name c c semi_elasticity if group == 1 & loss == 1 & post == 1, appendto("`tabname'.tex") rstyle(tabular)
    listtex elas_name c c taxterm_elasticity if group == 1 & loss == 1 & post == 1, appendto("`tabname'.tex") rstyle(tabular)
    qui listtex botrule if _n == 1, appendto("`tabname'.tex") rstyle(none)  
    qui listtex end if _n == 1, appendto("`tabname'.tex") rstyle(none)

end/*%>*/


**********************************************************************
* Table 9: Spikes in Capital Goods Shipments, Inventories, and Prices

capture program drop table_ppi /* %< */
program define table_ppi

	*gen industry specific timetrend to control for situations like computer related products. 
	egen mintime = min(time)
	gen timetrend = time - mintime
	gen yeartt = month - 1

	foreach x in 3 6 9 12 {
		gen d`x' = (month == `x')
	}
	gen dqend = (month == 3 | month == 6 | month == 9 | month == 12)
	egen industryid = group(industry)
	egen industryyear = group(industry year)
	eststo clear

	tsset industryid time
	gen lti_ratio = l.ti_ratio
	gen lti = l.ti_ratio
	replace vs = vs / 1000 
	
	eststo clear
	foreach y in vs vs_ratio lti lti_ratio ppi ppi_ratio {
		
		winsor2 `y', cut(1 99) replace
		eststo `y'1: qui: reghdfe `y' dqend /*c.timetrend#i.industryyear*/, ///
						  a(industryyear) cluster(industryid)
	    qui: estadd local indfyear "YES", replace  
		  
		eststo `y'2: qui: reghdfe `y' d3 d6 d9 d12 /*c.timetrend#i.industryyear*/, ///
						  a(industryyear) cluster(industryid)
		qui: estadd local indfyear "YES", replace  

		eststo `y'3: qui: reghdfe `y' dqend c.timetrend#i.industryyear, ///
						  a(industryyear) cluster(industryid)
		qui: estadd local indfyear "YES", replace  
		  
		eststo `y'4: qui: reghdfe `y' d3 d6 d9 d12 c.timetrend#i.industryyear, ///
						  a(industryyear) cluster(industryid)
		qui: estadd local indfyear "YES", replace  

		eststo `y'5: qui: reghdfe `y' dqend c.timetrend#i.industryid, ///
						  a(industryyear) cluster(industryid)
		qui: estadd local indfyear "YES", replace  
		  
		eststo `y'6: qui: reghdfe `y' d3 d6 d9 d12 c.timetrend#i.industryid, ///
						  a(industryyear) cluster(industryid)
		qui: estadd local indfyear "YES", replace  
		esttab `y'1 `y'2 `y'3 `y'4 `y'5 `y'6, keep(d*)     
	}  

	eststo m31: qui: reghdfe  ppi  vs  , a(industryyear) cluster(industryid)
	qui: estadd local indfyear "YES", replace   
	 
	eststo m32: qui: reghdfe  ppi_ratio  vs_ratio , ///
					 a(industryyear) cluster(industryid)
	qui: estadd local indfyear "YES", replace   
		 
	local tableconfig cells(b(star fmt(%9.3f)) se(par fmt(%9.3fc))) ///
		starlevels(* .10 ** .05 *** .01) noobs label ///
		collabels(none) noconstant sfmt(%12.0gc %12.2f) var(25) style(tex) nonotes 
	label var vs_ratio "Shipments (\%)"
	label var vs "Shipment Value"
	label var lti "Lagged Inventory)"
	label var ppi "PPI"
	label var dqend "Month 3/6/9/12"
	label var d3 "March"
	label var d6 "June"
	label var d9 "September"
	label var d12 "December"

	esttab  vs_ratio1 vs_ratio2 lti_ratio1 lti_ratio2 ppi_ratio1 ppi_ratio2 m32 ///
		    using "ppi.tex",  `tableconfig' replace  nomtitles  ///  
	mgroup("Shipment (\%)" "Lagged Inventory (\%)" "PPI (\%)" , ///
	pattern(1 0 1 0 1 0 0 )span prefix(\multicolumn{@span}{c}{) suffix(})) ///
	keep(d* vs*) scalar("N Observations" "r2 R-Squared"  ) 
	
end /* %> */

********************************************************************************
* Main
********************************************************************************

capture program drop main
program define main

	capture log c
	log using "$outputdir/draft_kinky.log", replace
    cd $outputdir
	
	global controls1 mkt2bk res2at capx2ppe sale4_3 l_at
	global controls2 mkt2bk res2at capx2ppe sale4_3 ebitda2at l_at
	global controls3 mkt2bk res2at ebitda2at l_at
	
	global graphconfig graphregion(color(white) lwidth(large)) bgcolor(white)
	
    *************************************************************
    * Graphs
    *************************************************************

	*** Figure 1 ***
	load_analysis_data capx_trim5
	fig_sawtooth_eps

	load_analysis_data quarterly_capx
	fig_quarterly_capx

	load_analysis_data fyrchange
	fig_fyrchange

	*** Figure 2 ***
	load_analysis_data international
	fig_sawtooth_intl

	*** Figure 3 *** 
    * fig_spike_bandwidth   //uses internal data 
	
	*** Figure 4 *** 
	load_analysis_data capx_trim5
	fig_uncertainty

	load_analysis_data capx_trim5
	fig_autocorrelation

	load_analysis_data responseq4
	fig_cumlevel_qu
	
    *** Figure 5 ***
   * fig_bunching    //use internal data 
	
	*** Figure 6:  ***
	load_analysis_data model
	fig_model_agg_spikes 
	
	load_analysis_data model_eventtime
	fig_model_cumlcompare

	load_analysis_data model
	fig_model_spikes_by_earn_bar, bins(10) 

	*** Figure 7 ***
	load_analysis_data m3
	fig_seasonal_agg

	load_analysis_data m3_ppi_nde
	fig_seasonal_agg_prices 

	*** Figure 8 ***
	load_analysis_data supplier
	fig_supplier 

	*** Figure 9 ***
	load_analysis_data lending
	fig_lending 

	load_analysis_data RateWatch
	fig_loanrate   

    *************************************************************
    * Tables
    *************************************************************

    *** Table 1 *** // no data used for this table

	*** Table 2 ***
	load_analysis_data capx_trim5
	table_sumstats

	load_analysis_data international
	table_sumstats_intl

	*** Table 3 *** // uses internal data
	
	*** Table 4 ***
	load_analysis_data capx_trim5
	table_constraints 
		
	*** Table 5 *** // uses internal data

	*** Table 6 ***
	load_analysis_data capx_trim5
	table_budget 

    *** Table 7 *** see table7_etrcalculations.xlsx in replication package
	
    *** Table 8 ***
    load_analysis_data model_did
    table_model_elasticity

	*** Table 9 ***
	load_analysis_data ppi_table
	table_ppi 

	cd $localcodedir
    log close
end
